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Abstract 

We present a perturbation theory by extending a prescription due to Feynman 
for computing the probabihty density function for the random flight motion. The 
method can be apphed to a wide variety of otherwise difficult circumstances. The 
series for the exact moments, if not the distribution itself, for many important 
cases can be summed for arbitrary times. As expected, the behavior at early time 
regime, for the sample processes considered, deviate significantly from diffusion 
theory; a fact with important consequences in various applications such as finan- 
cial physics. A half dozen sample problems are solved starting with one posed by 
Feynman who originally solved it only to first order. Another illustrative case is 
found to be a physically more plausible substitute for both the Uhlenbeck-Ornstein 
and the Wiener processes. The remaining cases we've solved are useful for appli- 
cations in regime-switching and isotropic scattering of light in turbid media. It is 
demonstrated that under isotropic random flight with invariant-speed, the descrip- 
tion of motion in higher dimensions is recursively related to either the one or two 
dimensional movement. Also for this motion, we show that the solution heretofore 
presumed correct in one dimension, applies only to the case we've named the "Shoot- 
ing Gallery" problem; a special case of the full problem. A new class of functions 
dubbed "Damped-exponential-integrals" are also identified. 



Key words: Random Flights, Brownian Motion, Feynman Diagrams, Perturbation 

theory, Path-Integrals, Stochastic Movement 

PACS: 02.50.Ey, 05.40.-a, 36.20.Ey, 42.68.Ay, 89.65.Gh 



Email address: tsh@MatheinaticusLabs.com (S. T. Hatamian). 
Preprint submitted to Elsevier Science 2 February 2008 



1 introduction 



The understanding of Brownian-motion has been of fundamental import on 
many levels in pure and applied physics. For example, the propagation of 
light in gaseous media has long been important in atmospheric and stellar 
physicsdll). In chemical-physics the topology of polymers(,2) is of great interest. 
Yet, unlike topics of comparable ubiquity, Brownian-motion has immediate 
applications in the human realm, ranging from the physics of finance(j3;0), to 
medical imaging(0), to name only a very few. 

The diversity of the applicants has resulted in a plethora of methodologies 
dehneating the evolution of this ubiquitous topic. Among these techniques, 
the path integral formulation of the stochastic movement, while often more 
imaginative, is fairly under-represented. Meanwhile, the standard methodolo- 
gies based on differential equations, often lead to difficult and unintuitive 
equations as soon as one gets away from the limiting situations such as in- 
finitesimal mean-free-paths, or long-time regimes. 

For example, the propagation of light through a stochastic medium is tradi- 
tionally described in the context of astrophysics by a Boltzmann transport 
equation for the distribution of intensity /(r, Q, t) in a heuristic Radiative 
transfer theory (0). However, since the general analytic solutions are unknown, 
one resorts to the diffusion approximation which can be shown to arise out of 
the Radiative transport equation in the limit of large length scales X»Lq, 
where Lq is the transport mean-free-path of light in the medium(j^; . 

In this same vein, there is considerable interest in the description of multiple 
light scattering at small length scales (X~Lo) and small time scales AT~1 
where 1/A is the transport mean free time, both from the point of fundamental 
physics (Isj) , and from the point of many applications such as medical imaging, 
etc. figtliol). It has been experimentally shown that the diffusion approximation 



fails to describe phenomena at time scales where AT<10 (jlll). Moreover the 
diffusion approximation, which is strictly a Wiener process, for the spatial 
coordinates of a particle is physically unrealistic. It holds only in the limit of 
the mean free path Lq— >0 and/or A — oo or the speed of propagation c— >^oo, 
while keeping the diffusion coefficient {D ~ c^/A) constant. 

It is somewhat surprising that developing better alternatives to the diffusion 
approximation has had to wait until the 1990's. For a particle moving with 
invariant speed in a one-dimensional disordered medium, it has been known 
for decades that the probability distribution function P(X, T) for the displace- 
ment satisfies the telegrapher equation exactly (0). However, generalizations 
of the telegrapher equation to higher dimensions ([lih have been shown not 



to yield better results than the diffusion approximation (|l5l: Indeed we 
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shall show (in section VI) that even the forementioned one-dimensional exact 
solution does not fully describe isotropic scattering, and additional diagrams 
are necessary to allow true isotropy. For a concise account of the recent work 
in light scattering work see 



In a similar mooring is the description of price movements in openly traded 
markets, where the diffusion approximation has been taken as gospel since 
the inception of theoretical finance. Meanwhile it has been known for over 
a centurv(|l7l) that the actual probability distributions of the market prices 
deviate significantly from those expected based on the diffusive (Wiener) pro- 
cesses. More directly, there has never been any concrete empirical reasons to 
believe a property such as A — >■ oo holds for the movement of openly traded 
market prices. With A~ for the typical gasses which concerned the 

developers of the early physics literature (P), the approximation was more than 
adequate. On the surface it appears that this notion was taken from physics 
without proper care to ask if the underlying assumptions are satisfied in the 
movement of market prices. We submit that addressing this "misapplication" 
is well overdue. The recent wave of exploratory work on the notion of "regime 
switching" (Q displays the burgeoning dissatisfaction with the gospel of 
diffusion in econophysics in general. 



Eleg ant prescriptions for the solution of the Fokker-Planck equation exist with 
(ji^ and without (21; 2^) resort to path- integral techniques which address the 



problem of classical stochastic probability density functions for arbitrary pa- 
rameters and time scales. However, techniques generally become cumbersome 
when applied to motion not prescribed by analytic equations. Alas, many 
problems, particularly outside the realm of physical systems, involve discrete 
and/or non- smooth motion which are easy enough to conceive but not easily 
expressed analytically. The method we outline is arguably better suited for 
this type of problem than other, more sophisticated, theories. 

In this report we present a technique based on a method invented by Feynman((2 
allowing notions of Feynman diagrams for the random-flight movement. The 



technique of path integration (j2J), also known as sum over histories, is partic- 



ularly attractive in the case of discrete stochastic motion because of the vivid 



imagery it evokes, literally illustrating the mathematical processfl25f). 



In section-2 we summarize the groundwork on which we will build our method. 
In section-3 we present the full solution to an illustrative problem set up by 
Fevnman(|23l) but only solved for small displacements. Section-4 illustrates 
how extraneous features such as noise or force fields can easily be integrated 
in the "propagator". As the final practice problem, section-5 presents the 
solution of an important variation on the problem of section-3. This variant 
can be considered a microscopically plausible alternative to the UhlenBeck- 
Ornstein process. The ensuing two sections, consider the problem of random 
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flights under invariant speed, where the majority of the results presented are 
novel. We then conclude with two specific applications of selected results in 
modelling flexible (polymer) chains and stock prices. 



2 Feynman's prescription for stochastic processes 



In this section we review the prescription(l23l) for computing the probability 
density function for a random process using functional integrals. Suppose that 
f{t) describes a realization of a random process as function of time. Called 
a "probability functional", P[f{t)], then gives the probability density of a 
given realization. If we think of f{t) as an ordered collection of discrete values 
/i, /2, /n in the limit of n^oo, the probability of finding a realization f{t) is 
given by P[/(t)]P/(t) = lim™ P(/i)P(/2)...P(/n)c?/irf/2..rf/. within a hyper- 
neighborhood given by the latter expression. The conjugate "characteristic 
functional" is computed by: 

_ /e-/^W/W^^P[/(t)]P/(t) 



Feynman(j23D showed that despite the challenging appearance of this expres- 
sion, it might be possible that under certain practical conditions the discretized 
form of P[f{t)]Vf{t) can yield to either term by term integration, or known 
forms of path-integrals. Suppose that each random event has a time signature 
given by g{t). Then the full realization of n-events over a fixed time interval 
T is simply f{t) = J2]=ig{t — tj). If the events are random and uniformly 
distributed in time, then the probability of each event occurring within the 
interval dt is dt/T. Inserting f{t) and P[/(t)] into Eq.(l) and noting the iden- 
tical form of each integral we get: 



^ n 

^ng[k{t)] = 1 exp j k{t)g{t - t,)dt 



\0 



J. . n 

gi/ k{t)g{t-T)dt'f[_ 

T 



dti 



dtn 

T 



(2) 



where we have used r as a stand in for any given tj. For typical applications, the 
shape or at-least the amplitude of each event is often not fixed. We therefore 
consider the case g{t) = au(t) where the shape u(t) of the event is fixed but 
its amplitude a varies according to some probability density p{a). We must 
now average ^ng[k{t)] over all values of a. Being that the latter is of form 
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{(f)g[k{t)])'^ . and that the events are independently distributed, we can carry 
out this average over all ^''s separately for each component (pg-. 



i<(pg >g-- 



J —J exp iaj k{t)u{t — T)dt 



p{a)da. 



(3) 



Observe that the inner expression is a fourier transform of p{a) where the 
frequency variable is the integral in the exponent. The fourier-transform being 
a characteristic-function, will be designated as W[J k{t)u{t — r)dt\. Finally, if 
the number of events over T is governed by a Poisson distribution, we can 
readily average over all possibilities of n: 

mt)] = j:r-Te-''. (4) 



Here n ~ XT is the average number of events expected over T. The sum is 
simply an exponential: 



= e-(^-^)" =exp 



1 

-xj(l-wj k{t)u(t - T)dt ) dT 



(5) 



An important special case of the process described by Eq.(5) is when the 
signal-event has an extremely short duration in time: u{t) = S{t). Then the 
characteristic functional becomes: 



= exp 



J 

-A J {1 - W[k{T)]) dr 



(6) 



3 The Gaussian Scattering Process 



Using the characteristic Eq.(6) we can compute various moments of the sig- 
nal profile f{t). However, the true utility of this approach is to discover the 
resulting distributions of dependent variables which are driven by the events 
in f{t). In this, and the following two sections, we will investigate cases of 
movement of free particles of unit mass, which undergo a prescribed change 
in their velocity for every event in f{t). Initially the particles have a delta- 
function distribution in both speed and spatial spread. In between events, the 
movement of the particle is governed by its appropriate equation of motion. 
Our objective is to find the resulting spatial distribution of particles after a 
given time T. 
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At every event the particle acquires a change in momentum which is a ran- 
dom selection from a gaussian distribution of spread s and zero mean. The 
characteristic function W{u!) of the gaussian random force- function profile is 
simply another gaussian. When inserted into Eq.(6) we have: 



mt)] 



—AT 

e exp 



(7) 



As discussed in the previous section, the probability density of the force func- 
tion P[f{t)] is given by the fourier transform of $[A;(i)], as given by Eq(l): 



In order to proceed, we discretize the time interval [0, T] into N, regular small 
intervals of duration 77. With the understanding that N^oo we can write: 



Pff,:! = e 



-AT 



exp 



Xrje 



-ivhfi 



(9) 



Finally, we expand the first exponential in a Taylor series and carry out the 
individual fourier transforms: 



N 



Xr) 



(10) 



whereupon wc recognize fiT] = li as the impulse imparted at time tj. The 
leading (zeroth-)order in is a product of (5-functions which we'll call a 5- 
functional 5[I{t)\. The value of the functional is zero unless the function I{t) 
is zero for all t. The zeroth-order term 

Pom] = 5[m] (11) 



corresponds to the ballistic path where the particle experiences no collisions 
and hence experiences no impulses. The 0{ri) terms are comprised of all 5's, 
except at time tj where a single term of 0{ri) contributes. There are N such 
terms (one for each tj) comprising a Riemann sum. In the limit of large N this 
sum becomes the integral: 

1 

Pm)] = -7^jd^ e-^W^/^^^5[/(t)]. (12) 
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The functional is defined as the (hniTv^oo) product of 6{I(ti)) for all tj 

except ti = t. This term describes the path where there is only a single collision 
at time t. It is easy to show that if the scattering profile had an offset, then it 
would change the exponent of the integrand to I(t) — Iq. Higher order terms 
in r] are easily obtained: 

Pn[m\ = n —/^ I dUe-'^'^^"''^'5[I{h), I{t^)] (13) 
,_i sy Z-K J 



where tn+i = T. Upon collecting the orders we get the expression for the 
probability density functional of Eq(lO): 

oo 

In order to obtain the probability density of a given output position P{X) 
we must first derive P[x(t)] from P[/(t)]. We then will sum over all paths 
x{t) which satisfy the boundary conditions: x{t = 0) = Xq ; x{t = T) = X; 
xlt = 0) = vo;x{t = T) = V. 

The final sum over paths can easily be done because of the insight we have 
acquired by recognizing the various orders in the series of Eq.(14) as paths 
with a given number of events. The zeroth order term is the easiest; it is 
simply the sum over all paths which experience no impulses, of which there 
is only one. However, to proceed more formally we will utilize the following 
obser^tionQ. As long ^ the relation between / and x is linear (e.g. I oc x), 
we can be sure that any Jacobian resulting from the change of variable: 

P^[I{t)]VI = PMt)]1^x (15) 

is a constant, and if skipped, affects only the normalization of the final answer. 
The normalization of the resultant distribution must be ensured regardless. 
Thus, from Eq.(13), and (15): 

Po{X,T) = J 6[I{t)]VI = J 6[x{t)]Vx = S{X- v^T). (16) 

The above relation (arguably the easiest path integral known) establishes the 
"propagator" for the process. Using Eqs.(13), (15), and two applications of 
Eq.(16), the first order term is the sum over all paths with a single event and 
is given by: 
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T +00 +00 

Pi(X,T) = —= f dti [ dxi f dh 

sy^'K J J J 

—00 —00 

X 5(X-((^;o + /i)(T-ti) + xi))e-^i/'^'5(xi-^;oti-xo). (17) 

This expression is easy enough to evaluate using the rule b{ax) — 5(x)/|a| and 
we will do so shortly. 

3.1 Feynman Rules 



Upon writing the sequence of terms Pi, P2, ■■■ we can see the that the unevalu- 
ated integrals (such as Eq.(17)) can be represented by diagrams. The diagrams 
comprise propagators, interaction points and the necessary factors to integrate 
over intermediate variables. In one-dimension the n^^ term is given by: 

Pn{X, T)= j d^Xn ... 1 d^XiKfnVn . . . K21V1KU (18) 



where / d^Xi = dxi J^^^ dti. The specifics of the problem are contained in 
the propagator K and the interaction profile V. Each term in the above equa- 
tion, and the final summation thereof can be represented diagrammatically as 
in figure- 1. 

For the case of a free particle undergoing a scattering process (which conserves 
momentum), the propagator is: 



1 / — * 

l^i+l '^i •' \ l^i+l k=l 



vA. (19) 



(Note: we will only treat the case of free particles here, so we will not clut- 
ter the super/subscript notation to that effect.) The interaction points Vf ~ 
/2s contain the gaussian profile of the i*^ event. 

The propagator can be generalized to other types of motion between collisions 
by simply replacing the classical equation of motion as the argument of the 
5-function or other appropriate kernel. 

At higher orders, the computation gets increasingly difficult and a recursion 
relation between orders of P„ does not exist in general. This is true in the 
present case of ^'s-process. For example to get Pf * we must repeat the integra- 
tion over X2 because K^^ involves X2. Therefore a more complicated integral 
over X2 must be performed before the integration over x^ is done. Despite this 



8 




Fig. 1. The diagrammatic representation of brownian motion through a medium. 

complication, the computation of P^*(X, T) can be carried out for arbitrary 
n, giving: 

Pf (X, T) = -1= fdt^ . . . / dt.i, exp (20) 
where the "interim variance" is given by: 

n 

ii'' = s'Y.(T-Uf. (21) 
j=l 



As we will see below, the summation over the orders can be easily carried out 
in the fourier domain. 

Should for some reason, the effect of scattering medium be artificially stopped 
after a fixed number of events, then the probability distribution is given by a 
new class of functions which we have named: "Damped exponential integrals" . 
These functions are described in the appendix-A. 



3. 2 Normalization 



We can directly demonstrate the normalization of P^'^{X,T) by integrating 
Eq.(20) over all X, and inserting the resulting terms in Eq.(14). Moreover, it 
is interesting to note the temporal population of each "generation" is: 



-XT\n 



iV.(r) = I P„(X, T)dX = ^^e-^^. (22) 



Reflecting our initial assumption, this is simply the Poisson distribution. Sum- 
ming over all n shows normalization. The population of generation n will rise 
and fall according to the above relation each with a peak at T = X/n (figure-2). 



Another way to verify the normalization of the gaussian-scatter solution is as 
follows. The quantities Nn{T) must satisfy the simple coupled rate equations: 
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Fig. 2. The relative population Eq.(22) of generations n, being the number of 
events they have experienced. Curves shown here are for A = 0.1. The peaks occur 
at T = A/n. 

A^o = -ATVo (23) 

It easily verified that the solution for these equations is given by Eq.(22) 

3.3 Moments and Kurtosis 

We will shortly demonstrate that the exact characteristic function for the 
^-s-process ^gs{k,T) can be obtained. Therefore all moments can be readily 
obtained by simple differentiation. However, the ability to obtain the charac- 
teristic function for a given process in closed form is not guaranteed. Hence 
we will demonstrate the direct computation of the moments by summing the 
moments of each order. We shall specialize to the case where all initial settings 
Xq,Vq,Io are zero, we compute The m*^ moment of the n*^ order distribution 
and then sum them according to Eq.(14). For the second moment this is: 

V y 

where Ti—{T — ti) . The discrete sum is nr^ by symmetry. The integral on 
is simply T^/S and the remaining n — 1 integrals of measure unity give T""^. 
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Inserting the result in Eq.(14) gives: 



<X'>,s= Is'T'e-'^ f = Is'XT^ (25) 



n=l 



{n-l)\ 3 



in agreement with the first order computation in (j2 

Based on the statement of the problem we can expect < V^>o(. s^T. Tradition- 
ally the computation of Brownian motion is within a system in equilibrium 
with a well defined temperature. In that case (jlj; l25i) the thermal equilibrium 
constrains the parameter s to the coefficient of friction and the temperature, 
and all three to <V^>. In the present case, no friction exists, therefore nor a 
finite temperature as evidenced by the ever increasing <V^>. The cubic de- 
pendence of the variance on observation time is reminiscent of the diffusion of 
tracer particles in turbulent flows. Richardson (EgI) was able to produce such 
a behavior for in the context of continuous diffusion only using a diffusion 
coefficient which depends on position as R^^^ (alternatively a time depen- 
dence of form could also do it(|27l)). More recently, it has been shown (|28l) 
that a Levy-distribution of waiting times can obviate the need for a space or 
time dependent diffusion coefficient in getting several families of supra-linear 
variances, one of which includes T^. The above result shows that a Poisson 
distribution can also give for the variance. 

A similar but more tedious computation for the fourth moment produces: 
X,s^<X'>,,= s'XT\^ + ^). (26) 



The first term is the result of cross-terms like (T — ti)'^{T — tj^ in the time 
integrals, whereas the second term is the result of direct-terms like(T — tj)^. 
Ostensibly, it is the presence of the interference terms that cause deviations 
from a normal distribution. Although as expected (from the central-limit- 
theorem) the deviation is transient. To quantify this, we compute the kurtosis 
for the process: 

""''^i^'^^^- ^^^^ 



Thus the distribution approaches normality in time. It is, moreover, worth 
noting that this result is independent of the parameter s. 
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Fig. 3. Demonstration of The Kurtosis (Fat-tails) as calculated for the ^s-process 
at AT = 10 versus the computer simulation of the ^s-process (The noisy line). 
The fit of the latter to a normal distribution is also shown. As expected from the 
central-limit theorem, the fit to a normal distribution progresses over time from the 
peak to the tails. Unlike the exact solution, the ballistic peak is not present in the 
normal distribution. 

3.4 The Characteristic Function 

Having arrived at a solution for P^^(X, T) (Eq(20)), the characteristic function 
is easily found by exploiting the symmetry in ^^"^ in the fourier transform 
of each P^^. The symmetry allows us to convert the sequence of connected 
integrals into the n*'* power of a single integral from to T with an added 
1/n! pre-factor. We can then sum them according to Eq.(14). The result is: 

M*.T)^.-exp(;fA£^;(^f|)). (28) 

Both the second and fourth moments can readily be verified by repeated dif- 
ferentiation of ^gs{k, T) with respect to k. The distribution Pgs{X, T) is shown 
in figures-3 and 4 versus a computer simulation. 

If it should happen that the deterministic equation of motion of the projectile 
is of the general type, x{t) oc fot"(0 < a < 1), then the result is a simple 
redefinition of Eq.(21) and for certain rational values of a the characteristic 
function can be easily computed. The result is that the variance will be pro- 
portional to T^+'^°' and fourth moment will be such that the kurtosis remains 
proportional to 1/T, independently of both a and s. This is consistent with 
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Fig. 4. Pgs{X,T) as calculated by a numerical simulation of the gs-process. The 
one-dimensional distribution is shown at different observation times T. For com- 
parison, at the corresponding time, the prediction of the theory (numeric inverse 
transform Eq.(28)) is superposed on the widest distribution. Only the amplitude 
of the theoretical curve has been adjusted to match the simulation. The spike at 
X = is the ballistic peak corresponding to the zeroth order theory (not shown). 

the Central Limit Theorem. 

3.5 Higher Dimensions 

For the case of a homogeneous and isotropic medium, and where the scattering 
profile has no preferred axis, it is easy to verify starting from the fundamentals 
of the theory, that the answer is found by replacing the scalar wave-number 
k, with the magnitude of the wave-vector. In case such rotational symmetry 
exists, by applying the inverse fourier transform we arrive at the probability 
density integrated over a D-dimensional spherical shell at R: 

oo 

^^^^ = (MpT^ I/^^^^ " ^' "^^^^ ^^^^ 

where, R = |X|, and J{D/2 — l,kR) is the Bessel function of the first kind 
and of order {D/2—1). Thus, if the g's-process is applicable, we can insert the 
D*^ power of the characteristic function directly into the above relation. 

Alas the g's-process is not applicable to photons because the postulate of rel- 
ativity is not satisfied by the infinite tail of the Gaussian interaction kernel 
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and the change in momentum after each scattering event is not independent 
of those along other dimensions. We will treat such constrained interactions 
separately. 



4 The Effect Of Noise 



If it should be that instead of simple traversal of the medium the particle 
diffuses in between collisions ("noisy paths"). We can show that insight pro- 
vided by the imagery of paths allows us to see our way through this added 
effect. Suppose Pgsd{y,T) is the probability density of arriving at {y,T) when 
diffusion is present, arriving at position y can proceed along an infinite num- 
ber of paths which would have culminated at different a;'s if diffusion were not 
present. If we let y = x+z then z is the component of motion due to noise. Thus 
the Pgsd{y, T) is the sum of all pairs (x, z) which satisfy the above constraint: 



Pgsdiy,T) = jdx Jdz Pgs{x,T)Pd{z,T)5{y - {x + z)) 

= JdxPg,ix,T)Pdiy~x,T). (30) 

From here it can be easily shown that the effect can wholly be incorporated 
into a new propagator K^'^. This is what we would have expected based on 
the previously alluded path imagery. The new propagator (itself obtainable 
by the path integral method (jl^) is simply a broadened form of 



where is the diffusion coefficient. The computation of Pgsdiy, T) proceeds as 
before. Computation will show that in presence of diffusive movement between 
scattering events, the diffusion acts in parallel to the process and thus simply 
adds a net term to the variance of the process. Thus, the "interim variance" 
is: 



■gsd? _ tgs2 



And we find for the full variance: 

l^gsd =< y"^ >gsd= ^S^AT^ + Cr^T; 



(32) 



(33) 
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the fourth moment: 

o o 



(34) 



and the kurtosis: 

^ Qs'XT 9 , . 



We observe that over fong times, the effect of diffusion (noise) becomes in- 
significant. Finaify the characteristic function is found to be: 

%sd{k, T) = e-(^-^^V2)T (^yi^ Erfl^_^^ ) . (36) 



In this section we have seen that diffusive intra-event movement can be repre- 
sented by a change in the functional form of the propagator. Other effects such 
as that of an external force field can be added to the theory by incorporating 
the equation of motion in the argument of the propagator. For an absorptive 
medium, a damping factor in the propagator will be necessary. 



5 The Gaussian Reset Process 



Consider if instead of scattering, we characterize the events as "resets" in the 
particle's momentum. It would be as if the previous momentum is lost at each 
event and a new one selected from a gaussian distribution of width s. This 
change implies only a different propagator: 

• = / S\ I - ^itiH^ ] (37) 



Utihzing this propagator in the Eq.(18), and other rules we find the probability 
density function P^^{X, T) to be the same form as the Eq.(20) if only we use 
a different "interim variance" function as given by: 

n 

iZ" = s'Eit^+l-t^r. (38) 
i=l 



Following the same steps as in the 5'5-process, we can show that in presence of 
diffusive movement between reset events the diffusion acts in parallel to the 
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Fig. 5. The comparison of the variance as calculated from a simulation of the 
(^r-process (dots), compared to Eq.(40). The upper solid line is that of 2s^T/A 
which is the progression of variance in (Wiener process) applications such as in 
Black-Scholes options-valuation. At early times, the linear approximation grossly 
overestimates the gr-variance. 

process and thus simply adds a net term to the "interim variance" : 

C = C + ^^T. (39) 

Unfortunately, for this process, ^^^^ does not posses the symmetry which al- 
lowed us to compute the characteristic-function in closed form for the gs- 
process. Thus we must compute the various moments by summing the mo- 
ments of each P^'^{X, T) according to Eq.(14). This task, while a little tedious 
for higher moments, is straightforward. The variance is found to be: 

V,, ^< X' >,,= ^ (2(e-^^ - 1) + AT(e-^^ + 1)) . (40) 

In very early times the variance increases as which matches Vgg. Later, 

the momentum non-conservation in this process manifests its different char- 
acter resulting in slower spreading of the distribution. For this process the 
variance becomes linear in T in the long time limit (AT ^ 1). 

The fourth moment is: 

"X.gr =^ ^ ^ gf 
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^[e-^^{12 + 18AT + QA^T^ + 2X^T^) - 12 - 6AT + SA^T^]. (41) 
A 



At early times, the fourth moment also behaves like Xgs increasing propor- 
tionally to T^, but later it settles into a parabolic increase. The kurtosis for 
the gr-process has a steep fall off at early times, but then it too settles into 
~ 1/AT descent: 



AT -4 

Tlarge (AT — 2) 



Kgr Pf^ 2 ^^^ (42) 



Observe that the eventual linearity of the variance is reminiscent of the Uhlenbeck- 



Ornstein ( [/0-)process (see e.g.(l2ll).ch.3). The t/O-process is random flights 
(i.e. the Wiener-process) in the presence of continuous friction. However, con- 
tinuous friction may not be a plausible mechanism at the "atomic" scale, 
where by an "atom" we may mean: an atom of matter, a node on a network, 
an individual trader, etc. Hence it is tempting to to think of the gr-process 
as the microscopic alternative to the mesoscopic [/0-process. As a reminder, 
both the UO and gr-processes approach the limit of the Wiener process, if we 
turn off friction in the former (7 but 7^ finite - as defined below), or allow 
very large number of collisions per unit time (AT ^ 1) for the latter. 

For the [/0-process, in time, the variance approaches uuo = 26T/m'-f, where 6 
is the absolute temperature and 7 is the coefficient of friction. Thus to match 
this with the gr-process we must have: O/rwy = s^/X. if we arbitrarily set 
6/m = s^, and 7 = A, then a graph of looks quite a bit like that of (fig- 
5). However, it is not possible to exactly match the variance of the f/O-process 
to that of the gr-process for all times for any conceivable relation between 
(s,A) and (^,7). This stems from the fact that these parameters represent 
quite different physical meanings. Still it is interesting to note that (unlike the 
5's-process) the gr-process has an effective temperature due to the stabilizing 
effect of the reset mechanism which acts clS clS db form of friction. Thus the 
choice between the U O vs. gr-processes as the appropriate physical model for 
a given system lies with the end-user. 

We conclude by stating the moral of this section. Even if the probability 
distribution (or its characteristic function) may not be computable in closed 
form, in many cases, our method allows access to the exact moments of the 
distribution. 
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6 The Shooting Gallery 



To further demonstrate the versatility of this method we will take on a prob- 
lem where a solution by the "traditional" approach of stochastic differential 
equations is difficult. Consider the carnival game of "Shooting Gallery". A 
"target-duck" starts moving at the center of a plank with a fixed speed c to 
the right. The player is required to shoot the duck. On every hit the duck 
reverses direction and continues to retrace its path with speed c. After a time 
T, we need to know the probability distribution of finding the duck at a given 
position X within dX. As before, the number of shots is given by the Poisson 
distribution of mean AT, but otherwise uniformly distributed over [0,r]. 

The scattering profile (i.e. the distribution from which the next momentum 
change is selected) is now comprised of two 5-functions, each offset by the 
amount ±c, respectively. Each lobe of the scattering profile is used in alternate 
order. The alternation requirement is an instance of a situation where the 
description of the motion does not lend itself well to analytic description. 
Because there is a binary alternation between left-right symmetric profiles we 
will call this process the "symmetric-binary-delta-scattering" , or s&c?s-process. 

The propagator for this problem is the same as that of gs-process (Eq.(19)). 
Furthermore the result for the Pf^ when vq (Eq(17)) apphes if we replace 
the interaction kernel (gaussian profile) with S{Ii+2c). Higher orders can easily 
be written using the Feynman-rules for this case. Integration over intermediate 
positions can proceed as before and the analog to of Eq(20) is found to be: 

psMs^X,T)^ jdt,,... I X-(^;o-c)T-c^(-l)-(i,+i-t,)j(43) 



where to = 0, and tn+i = T. We consider only the initial condition Vq = c 
which results in cancellation the second term. Although higher orders of the 
above integral require a good bit of bookkeeping, the actual integration are 
straightforward if aided by computer. By extensive use of the rule, J^dr5(j- — 
t)f{T) = f{t)[e{b- t) - e(a- t)], we find Pf^'{X,T) to be: 



5{X- cT) 



n = 0; 



^ 2 B{X,cT) n = l,3,5,...; 



2 c" (n - 1) 

(x + cr)(c^r^-x^)^ 

2c"n!! (n-2)!! 



B{X,cT) n = 2,4,6,...; 



(44) 
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Fig. 6. The time evolution of the sbds-process for c = 1, A = 1. Here The zeroth 
order {6(X — cT) in Eq.(45)) balhstic peaks have been artificiany broadened for 
plotting. 

where n\\ = {n/2)\ 2"^^ (ra-even, e.g. 6!! = 6-4-2). The "light cone" is maintained 
by the "box-car" function: B{X, cT) = Q{X + cT) - Q{X - cT). Finally, the 
observed probability density is found by summing P^"^^ according to Eq.(14). 
To aid the summation process we can shift the dummy index n 2r;, + 2 for 
even terms and n — > 2n + 1 for the odds. The summed result is: 



P,,,,(X,T) = e-^^|5(X-cT) + 



B{X,cT) 
2c/A 



(X + cT) 
pc/X 



/(l,p)+/(0,p) 



(45) 



where I{k,p) is the modified Bessel function of the first kind and of order k. 
p is a measure of position within the light-cone, given by: 



We note that this is in agreement with (|l2f l derived by a laborious solution 
to a coupled set of "telegrapher's equations" . Figure-6 shows samples of Psbds 
over a relatively early time range for A = 1 and c = 1. 

The solution above and that in the figure pertain to the initial condition: 
P(vo = c)=l. In the event P{vo = —c) 7^0 then the solution can easily be con- 
structed by the superposition P {vo = c) Psbds {X,T) + P{vo = ~c)Psbds{—X,T). 

Finally, we could have modelled the problem with a reset type propagator as 
in Eq(37). It is easy to show that the resulting sMr-process will give the same 
answer as long as the initial condition vq = ±c applies, but not otherwise. 
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7 Isotropic Random Flight in D-dimensions at an Inveiriant Speed 



Extending the previous section's model to higher dimensions is of great prac- 
tical interest in physical systems. Our method facilitates the setup without 
difficulty. We shall denote the unit vector along the d*^ component of motion 
after the j^^ scattering event as ejd- As before, we assume that the particle 
moves at a constant speed c, and use the propagator of Eq.(19) for a free, 
momentum-conserving particle. After the easy integration over intermediate 
positions we obtain for the end-to-end propagator of the o?*'* component of the 
position: 



where to — 0, and t^+i — T. The above integrand is quite plain in its statement 
about the free particle, and could have been written without the need for 
integration over many intermediate coordinates. We can consider two possible 
initial conditions: 1. An incident beam where = c, or 2. A source emitter 
where — 0. In the former case we can select Vq = ceo resulting in cancelling 
the middle term, but we must maintain eo = (0, 0, Alternatively, we 

can shift the sum to start from 1 but remember that the position vector is 
given by R = {x,y,...,z — cti). The latter case is simpler because we have 
eo = which simply eliminates the middle term, the sum remains as is, and 
R— {x,y, z). Either way, the computations do not depend on the absolute 
value of the index j, hence this distinction between initial conditions does not 
come into play until the very end whence integrating over interaction times. 

According to the Feynman rules, we must now add the interaction kernels 
and sum over all allowable intermediate directions Cjd as well as times tj. The 
allowable states for the particle include only those which maintain a constant 
speed c. Therefore, for movement in Euclidean space, we must implement the 
constraint J2d=i ^jd — ^ ^'^^ J- Thus the full diagram for the N^'^ order for 
the D-dimensional scattering (Z)ds-)process is: 



N +1 D 



oo 



X V{e) 5(^d-cj:^^e,d{t,+,-t,)^5(l-Y.e^^^ (48) 

As a consequence of the nonlinearity of the unit vector constraints, additional 
normalization factors will be needed depending on the explicit dimensionality. 
Here we will consider isotropic scattering and thus set V{e) — 1. 
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In performing these summations we employ a technique which could also ap- 
ply to most of the problems we considered previously, but did not for the 
sake of illustration of the possibility of direct integration. The 5-function can 
be represented as an unweighted superposition of plane waves. This action, 
amounting to a fouricr transformation of the integrand, results in the compu- 
tation of the characteristic function instead the probability density. However, 
in this case the decomposition is especially necessary. Once decomposed, terms 
involving a single eja are collected and integrated separately. There is, how- 
ever, one further complication: An integration over ejd diverges unless the 
wave number corresponding to the second 5-function is always positive. For 
this reason we choose an alternate (if obscure) plane wave superposition given 
by, 5{x) = Re{J^ e^^^dK,/7T)—i/T:x, where the second term is understood as 
the principle value. Together with the ordinary decomposition of the first 5- 
f unction, the integrals over involving terms like i/irx will all vanish. The 
remainder of terms with Cjd each result in ^J-K/iKi exp (iA;Jc^(ij+i — tjf /Anj) 
where the nested products over j and d apply. In order to carry out the inte- 
gration over the Kj we perform the product over j. This conveniently results 
in the square of the magnitude of the vector k in the exponent. Thus we 
arrive at: 



1 /■ e 



where, Sj = kc(tj+i — tj)/2 and k = |k|. Applying the product over j to 
Re{(f)u), and the integration over all interaction times tj, we arrive at the 
characteristic function for the isotropic scattering after A^-events: 

^^'"(^)=n / Re{Mk:tj))dtj. (50) 

J=0 



Before proceeding, we note an important property oi (f)D'- 

'^^-^^ ^ -Ys^- ^^^^ 



That is, the respective characteristic functions for all odd and separately, all 
even dimensions are recursively related. For this reason we need to evaluate 
the expression in Eq.(49) for D—1, and 2 only, viz., 

0i(«)=e^-; (52) 
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and after normalization (i.e. 0(0) — 1), 



(53) 



where K[0, x) is the zeroth-order modified Bessel function of the second kind. 

The probabihty density is found by performing a D-dimensional inverse fourier 
transform on However, because the $^**'(A;) possesses rotational sym- 

metry, the angular integrals can be performed and the fourier inversion reduces 
to Eq.(29). 



7.1 One-dimensional Isotropic Scattering 



Although the motivation for computing the Ddis-process was the generaliza- 
tion of the shooting-gallery problem to higher dimensions, it turns out that the 
latter is not the same as the onc-dimcnsional (idis-)proccss. This can be most 
readily seen via the characteristic function (inserting Eq.(52) into Eq.(50)): 



^'^^^(fc, T) = n / cos{kc{t,+^ - t,)) dt,. (54) 
^■=0 



By setting = 1 in the above and fourier transforming the first order function, 
we can readily sec that the 5-function found after the integration over the 
interim coordinates in the sMs-process (Eq.43) is only one of four we obtain 
here. We further note that one of the extra terms corresponds to a first- 
order process where no velocity flip takes place at the event time ti. The 
remaining two are the parity (^— > —X) conjugates of the latter two terms. 
These features reflect exactly the characteristics which produced Eq.(48). That 
is we only required that the magnitude of the speed to remain constant at all 
times, but allowed all possible directions. We also allowed isotropic initial 
conditions which added the parity conjugate terms. Schematically, the Idis- 
process encompasses all diagrams in figure-7 (plus their parity conjugates), 
whereas at each order, the sMs-process includes only the left most diagram 
in the figure. In the Idis-process, it is as if after each hit, the duck in the 
shooting-gallery fiips a coin to decide whether to reverse direction or not. 

While there are 2^ diagrams at each order, only (N-l-1) are distinct. Other 
than the ballistic term, present at every order, we label the remaining N by 
the index j below. Explicit computation of the this sub-collection yields: 

pM^s(yrJ.. MrT ^^ , B{X , cT)m (cT + Xy (cT - X)^'^-' 
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Fig. 7. The diagrams that are included in the Idis-piocess. The mirror images must 
also be included, but are not shown 
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Fig. 8. Pidis at AT =1, 3, and 5. The ballistic peaks have been artificially broadened 
for plotting. 

B{X, cT) is the "Box-car" function defined previously, which maintains the 
light-cone. We can now insert this result in the summation formula of Eq.(14), 
and add the parity conjugates, resulting in: 



pidis(^^^ T) = -{e-^^/^S{cT - X) 



+ e-'^B(X,cT)Y: 

N=0 

+ X}. 



~ ^ A\^+' (cT-X)^ 



4c, 



F([-7V-l,-Ar],[l],p)} 

(56) 



Here p 



^z^i and [7], p) is the hypergeometric function of the 

continuous variable p. Sample behavior of this probability density function is 
shown in figure-8. 



If accuracy is not detrimental, the following simpler expression for the par- 
tial probabilities has the same second moment, and approximates the fourth 
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moment of Pjf^^{X,T) to within 5%, 



Plf'\X,T) ^ N\{N+2) y~T^) ^^^^ 



applicable for N > 1. For T < 5, the full solution (Eq.(56)) without the 
ballistic term and the light cone, may also be closely approximated (in units 
A = 1, and c = 1 ) with the form: const/T^ ■ exp{(T+ 3)(1 -XV(T+ 3)^)^/^}, 
where const ~ .0022. The approximation can be very good if one is interested 
only in small variations in T. There, the appropriate value of the const can 
be found by mere eyeball fitting. For larger T, a gaussian starts to become 
viable. 



7.2 Two- dimensional Isotropic Scattering 



In two dimensions, the characteristic function is found by inserting Eq.(53) 
into Eq.(50): 



N *J +1 



^^t{k, T) = t[ J J{0, kc{t,+, - t,)) dtj (5J 



i=o 



where J(0, x) is the zeroth order Bessel function of the first kind. The fourier 
inversion is made possible via Eq.(29) with the kernel /cJ(0, kR)/2Ti. For the 
source emitter configuration, explicit computation yields: 



i-r(i^,T) = ^^;,^^_;^, (c-T--i?-)- (59) 



applicable for > 1. Summation of orders via the summation rule Eq.(14) 
yields the full probability density. 



Note here, that the diffusion limit is found when cT ^ R rather than the 
usually assumed limit of AT 3> 1. The former, however, is in line with the 
implications of the central limit theorem, whereas the latter is more of a "rule 
of thumb", which must be used with some care. The above two-dimensional 
results agree with those given in (jl4^ . 
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1.3 Three-dimensional Isotropic Scattering 



The 3-dimensional process is likely of greatest practical interest. Below we 
specialize to the case of the source emitter. Using Eq.(51) we find i?e(03) = 
sin{sj)/sj. Inserting this in Eq.(50) and then into Eq.(29) results in: 



-pZdis 



{R,T) 



1 



dti °f sin{kR)sin{ckTj 



27r^Rc 



n/^/- 



(ck) 



N 



dk 



(61) 



where Tj = {tj ^i — tj). As expected, the zeroth order integral produces the 
ballistic peak: p^&s ^ §(^cT ± R). The first order integrals can be found by 
convolution of 03 with p^djs direct integration to give: 



Pi 



3dis 



e{cT-R) ^JcT + R 



cT-R 



(62) 



The result of the sine-transform in Eq.(61), proves inaccessible for > 1. 
However, the moments of the distribution are easily calculated to arbitrary 
order. While all odd moments are zero, the even moments are found by the 
application of (V^)"^/^ to for m = 0,2,4, .... For m = normalization 

can be verified as M^'^{Q) = r^/(47rA^!). For m > 2 The m*'* moment is: 



rp(N+m) 2"'/2 I^J\J + 1) 

47r2c"^ {N + my. (m - 1)!! 



m/2-l 
j=0 



(63) 



where (m - 1)\\=2'^/^T{^)/ ^ for even-m (e.g. 5!! = 5 ■ 3 • 1). Only the set 
ao{m) = m!^/(2'^(m/2)!) yields to analytic description. The coefficients aj{m) 
for up to the eighth moment are listed in the table below. 



aj{m) 






i=2 


J=3 


m=2 


1 








m=A 


18 


5 






m—Q 


1350 


1715/3 


175/3 




m—8 


264600 


137018 


22785 


1225 
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For iV> 1, Paasschens([l3) has proposed the expression: 



Wr).^£i|±il^fi-i^f-' ,04) 



for the probabihty density of the 3c?zs-process. This expression produces the 
required second moment exactly, the fourth to within 0.5%, and the sixth 
moment is approximated to within 1.5%. Hence an excellent approximation for 
many practical purposes. The total probability density function (via Eq.(14)) 



has also been approximated by (jlj): 



—AT 

P'^'^R) ^ ^6{cT - R) 



+ e-^^e(cT - ^) \4,,24a)3/^ e^3y 1 + 2.026/P3. (65) 
where, p^d = AT(1 - R^/^T^f'\ 



8 Sample Applications 



The transmission of photons in turbid media is of interest in medical imaging. 
During the 1990 's increasingly more successful attempts have been made to 
model the stochastic movement of the photons in turbid media (flil : Hot Is^ 



31l : l32l : |33|). Some of these have involved forms of path-integration methods 
while others, not. However, all of these reports have been limited by varying 
forms of approximation, limited dimensionality, and the like. Some of these 
approximations pertain to truncated orders of computation or other more 
subtle ones such as maintaining the photons' light-cone (js^) only on average. 
Nevertheless, for practical purposes, computations of highly forward-scattering 
seem suitable for applications involving biological tissues. As such, the findings 
are in reasonab ly g ood agreement within the precision of measurements as 



reported in (|30l : l34f ) 



The characteristics of these works have been recounted in a chronological 
narrative in (jlfih . For the case of isotropic scattering, graphical comparisons 
of several of these works to certain exact results can be found in fig-3.3 of (14). 

Here we will not consider the mathematically intricate anisotropic scattering 
application as it does not make a good illustrative case. The large body of 
literature in that realm, nevertheless, is good evidence of the struggles of 
usual approximations with early times (AT < 10) in stochastic motion. We 
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will however, consider two simpler popular applications of our results: polymer 
chains, and stock option valuation. 

8.1 Flexible Polymer Chains 

A minimal model of a flexible polymer is a chain of links of constant length L 
and total length x. By "flexible" we mean that each bond is free to assume any 
orientation in space as long as it remains linked to its two neighbors. Therefore 
the probability density distribution of such a chain is a special case of what we 
have already considered in the Ddis-process where c(tj+i — tj) =L = cT / (n + 1), 
such that n, the order of computation, is the number of joints. This also means 
that the often-difficult time-ordered integration becomes unnecessary. That 
is, back in Eq. (2) the probability of an event at tj is not uniform over T but 
restricted to the single instant tj = T/j. Thus we must use 6{tj — ■^^)dtj 

in place of leading to the removal of the time integration in expressions 
such as Eq.(50) for fixed-sized steps, and resulting in factors which only affect 
the normalization. Moreover, because of the discovery of the recursion rule in 
Eq.(51) we need only work out the characteristic function for only one and 
two dimensions. All higher dimensions can then be derived from these, using 
the recursion relation. 

The one- dimensional case is likely of limited interest but we provide the results 
here for reference. Using Eq.(52) and Eq. (50) (without the time integration), 
the characteristic function for a chain of I links is the power of cos{kL), 
and the effective fourier kernel is cos{kx) for kE (— oo, oo) and, xG {—oo, oo). 
Writing the cosine as exponentials, the transform is easily found as a collection 
of 5-functions weighted by binomial coefficients: 

pin^) = 

j,jz[4^)j{5{x-L{23-l))+5{x + L{23-l)))- £ = 1,3,... (66) 

l{,%y{x) + j;ty-?j(^^^^ £ = 2,4,... . 

It should be apparent that for large ^ the binomial coefficients will tend toward 
a Normal distribution enveloping the discrete spikes (Fig-9). In that limit, the 
spikes become indistinguishable and one can approximate P^'^^'^^x = mL) as a 
single term given by e~™^/^"'/-\/27m(j3^. 

The two-dimensional problem has a famous historvllssh but it is easily ad- 
dressed using our recursion relation. In two- dimensions, using Eq.(53 and 50) 
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Fig. 9. The probability density distribution of Eq.(66)) for £ = 21 links and L = 0.1. 
The peaks have been artificially broadened out for plotting. 

the characteristic function for a chain of i hnks is the i^^ power of J(0, kL), 
and the effective fourier kernel is kJ{0, kr) for both k and, r G [0, oo). Despite 
the fact that the same transformation in conjunction with the time-ordered 
integration could be worked out in the case of 2dis process, we cannot analyt- 
ically obtain the probability density function for the end-to-end distance for 
higher than second order (£>3). The zeroth order contains only a single link 
and obviously corresponds to r)/27rr, and the first order is found to be, 
Q{2L — r)/7c'^r^/4L'^—r'^. The second order (3-link) chain has the probability 
density for the end-to-end distance, K{A/y/L^) /n'^y/lJr, for r G [L, 3L], and, 
'K(^/ L^r / A) / ir'^ A, for r G [0, L]. Here, K is the complete elliptic integral of the 

first kind, and A = {L+r)^L + r)(3L — r)/4 is the same as the area subtended 
by the quadrilateral formed by the the links in the chain and the end-to-end 
vector of length r. As alluded, higher order distributions remain analytically 
inaccessible. But since we have the exact characteristic function, computing 
the exact moments at any order is straightforward. Below we provide up to 
the eight moment, for chains of up to five links. 



M/L™ 


m=2 


m=4 


m=6 


m=8 


n=0 


1 


1 


1 


1 


n=l 


2 


6 


20 


70 


n=2 


3 


15 


93 


639 


n=3 


4 


28 


256 


2716 


n=4 


5 


45 


545 


7885 



For large-£(~ n), the distribution has been found by a number of methods 
over the last centurvljssh to approach: 2r/(L^n)e~''^''^^". 
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The most useful case clearly being that of 3-dimension has been worked out 
by Kleinert . We find agreement with this result as follows. Specializing 
Eq.(61) to the problem at hand, we find the characteristic function for a chain 
of ^ links is the power of sin{kL)/kL, and the effective fourier kernel is 
ksin{kr) /2T['^r for both k and, rG [0,oo). Because the characteristic function 
is even in k we can convert the kernel to an exponential and then extend the 
lower limit to —oo. If we write sin^{kL) as combination of exponentials we 
find an integrand of form e*^^//c^~^. The full integral can then performed using 
contour integration resulting in(j3): 



(67) 



Here the large-£(~ n) limit can be shown to be: ^^\f^e ^r-^/aL^n^ 



8.2 Financial Options Valuation 



The problems considered here are especially relevant to the price movements of 
financial assets. Economists traditionally assume that the logarithm of prices 
have a normal distribution (js^), yet it has always been known this is only true 
for very large A. While this assumption is fine for the intended typical gaseous 
medium where A~10^''/sec, the equivalent reset rate for a liquid market is of 
order inverse days (j37l ). Figure-10, showing the observation-time dependence 
of the kurtosis over 32 years of S&P log-returns, clearly indicates a finite A 
and hence the inappropriateness of the diffusion approximation. 

The fact that one cannot fit a single A for different observation-intervals (T) 
(ranging from minutes to many weeks), suggests that more than one process 
govern different time scales(38), a notion so common, it is taken for granted 
by traders. With more than one Pgr curve one can arguably reproduce the 
necessary fat-tails as well as Levy-distributions (figure-11), if only as another 
alternative. 

In order to see some of the implications of finite A, we present briefly, the valu- 
ation of "options" under the gr-process. Option pricing using path-integrals is 
quite commonplace; seeljsil: liot IZH to mention only a very few. The com- 
putation of the worth of an (European-type) call-option for a non-dividend 
paying underlying asset was first computed by Black and Scholes (0). A call- 
option is a contractual right purchased for an agreed upon premium C, which 
gives the buyer a right to purchase a fixed amount of some asset at a later 
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Fig. 10. Comparison of the kurtosis for the daily (Jan-1970 to March-2002) 
log-returns of the S&:P500 index vs. a fit of our typical k ~ 2/ XT. For 
1/A = IS.Jdays the agreement is good for periods less than 30 days. 




Fig. 11. The distribution of log-returns daily (Jan-1970 to March-2002) S&:P500 
index vs. Pgr, and Levy distributions. While the Levy curve fits the fat tails 
much better, a superposition of Pgr^s can also fit the fat-tails. Conversely, the 
Levy-distribution does not have the ballistic peak; a feature extant in observa- 
tion though conveniently ignored by many investigators. Clearly, the correct phys- 
ical model which gives rise to the Levy-distribution would produce the ballistic 
peak(!:44). The correct relative amplitude of the ballistic peak to the smooth curve 
in Pgr has not been force- fit. 

time T from the seller at an agreed upon ("strike") price 5*. This is regardless 
of what the prevailing market value of the asset might be at T. 

Black and Scholes(BS) first assumed that the logarithms of the asset's price 
follow pure diffusion, with an undetermined diffusion coefficient o"^ (a.k.a., 
volatility). The distribution 6{lnP — IuPq) reflects the full certainty that the 
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Fig. 12. The progression of price distribution according to Black and Scholes 

current price is Pq- Thus the distribution will spread out in time into a Gaus- 
sian of variance a^T. They also assumed the log-price to drift linearly (but 
no faster) in time with an undetermined rate [i (figure-12). In this way they 
thought they could fix the price probability density projected for the time T 
into the future. 

The payoff of an option for the buyer would come about if Pt (the market price 
at time T), is higher than the agreed strike 5; the option expires worthless 
otherwise. If Pt > S, the buyer could exercise the said right and immediately 
sell the asset on the open market for a profit of Pt — 5*. Thus a fair price to 
ask of the buyer is the weighted average of all possible payoffs: 



—rT 



oo 



Cbs = -^^^ I (e^' " S)e (6J 



InS 



where r is the prevailing interest rate. The reason for the pre-factor e~ is 
that the buyer loses money by missing out on a steady interest payment since 
s/he has given up the cash to buy the call-option. BS then argued: if it is true 
that no-one can win consistently at speculating, the drift rate n must exactly 
be equalSiJ to the guaranteed inteiest rate r. I,r this way one of the two 
arbitrarily introduced parameters was eliminated. The integral can be written 
with limits xE [— oo, oo], if we write the integrand with a step- function. Eq.(68) 
can then be integrated by parts resulting in an Erfc (called in finance texts). 
The final expression is known as the Black-Scholes formula. 

Ostensibly, there are assumptions in the BS hypothesis which are not sup- 



31 



ported by observation. Conversely, there are many observations which are not 
incorporated in the BS model, such as Levy-like fat tails in the probability 



distribution(j^; |46|; |47l). The resulting discrepancies have become deep puz- 
zles in the realm of finance (jiih. One such puzzle manifested in the pricing of 
options is the phenomenon of "price- skew" . In options trading practice the 
price of an option cannot really be calculated using the BS formula, even after 
eliminating fi, since the "volatility" (a^) is not known. This parameter may be 
inferred from the variance of the detrended price series of the underlying but 
using it does not give option prices that match the observation. Conversely, 
if the observed option values are inserted into the BS formula, the implied 
volatility cr^ , resulting from solving: 

Cobs = CBs{S,T,Po,r-a), (69) 



is not the same as that inferred from observation. In fact, using the observa- 
tions from different option series (expirations T or strike prices S) one gets 
different values for cr, whereas Cbs only allows a unique value for all options. 
Ostensibly, there are many more variables which affect Cobs in the real world. 
The existence of any additional variables or control parameters immediately 
implies that for any given Cobs, the solution for a is no longer unique. If we 
take the gr-process as the microscopic basis for the oft-presumed continuous 
price diffusion (Wiener) process, it provides just such a parameter in the form 
of A. 

It will be illustrative to compute the option skew implied by the gr-process. 
We will do so only in an approximate way so that we need not compute the 
full Pgr- This will also clearly demonstrate that the "fat-tail" of the distri- 
bution resulting from finite values of A is directly linked to the option-skew 
phenomenon. We will emulate Pgr by matching its kurtosis to that of a super- 
position of two gaussians (gausslets): 

2 2 

_ X _ X 

p 2!^^ p 2/^2 

We can now solve for {1^1,2/2} by setting equal the variance and the fourth 
moment of G2 to that of Pgr- For illustration sake we can take the intermediate 
time (AT > 10) limit of Ugr Eq. (40), and Xgr Eq- (41). In this approximation 
the solutions are: 

1/1,2 = ^ (2AT - 4 T 2v'(2Ar-8)) . (71) 

We can see that if we attempt (as traders do in relying on the BS-model) to fit a 
single gaussian to the observation, we are likely to be accepting only one of the 
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Fig. 13. The implied volatility (ordinate) vs. Strike price (abscissa). The 
two-gausslet (G2) representation of Pgr is used in place of Cobs in Eq(69). Parameter 
settings are: AT = 10, Pq = 50, s = 0.1, /i = 0, r = 0.05. While the features of skew 
plot remain the same, their location and span vary to a large extent for different 
parameters. For this reason one could get a smile, frown, or smirk relative to the 
at-the- money point S = Pq. The easiest way to obtain a smirk or a frown is to allow 
/i 7^ 0, this is consistent with the findings in (j49h. 



above variances. It is manifest that the implied volatility has a dependence on 
the time to expiration T other than the traditional linear term. To demonstrate 
the same for the strike price S we use the full form of Ugr and Xgr and solve 
for the implied gausslet variances {i'i,i'2}- This refined procedure allows us 
to get to times as low as AT ~ 2, but to go further we must match more 
moments. We produce an option valuation formula as the half the sum of two 
Black-Scholes type expressions but using {vi,U2\ as the respective variances. 
Settings this approximation to Cbs and solving for a for different strikes gives 
us the the skew curve as shown in the figure-13. 

The application of the method of path-integrals to financial derivatives' valua- 
tions is a very natural approach. Unlike the case of European (path-independent) 
options treated above, many options have values which are path-dependent. 
Hence the valuation must take place as part of the path-integration itself 
(|3i;H|4i|). 
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A Damped Exponential Integrals 



If it should happen that the gs (dj-piocess is terminated after a fixed number 
of events n, then the distribution (for the truncated-gaussian-scattering with 
diffusion) is given by a set of integral functions: 



pt9sd(^^^ T) = Dein{x, T; A, s, a) 



-\tr. 



dt 



Jdtn-i... Jdti -^exp - 

\ 



X 
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(A.l) 



where Dei represents the name "Damped Exponential Integral" . The "interim 
variance" is defined as before in Eq.(21). The normahzation factor is given by: 



Dei 



A" 



(A.2) 



Thus DetQ is the familiar normalized Gaussian of width V cr^T. Further for 
example, Deii is a superposition of Gaussians of all variances from cr^T to 
(a2 + s^)T: 



Deii{x, T; A, s, a) 



^2n{l - e-A^) J 



I 



dt 



-At 



X 



exp 



X 



^^s^^T-ty + a^T) V 2s^{T-tf + a^T^ 
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